// eee.Sheffield.PZ.Imaging
//
// Copyright ?Ping Zou, 2006
// sg71.cherub@gmail.com

using System;
using System.Collections.Generic;
using System.Text;
using eee.Sheffield.PZ.Math;
using System.Drawing;
using System.Drawing.Imaging;


namespace eee.Sheffield.PZ.Imaging.ImageFilter
{
    public class FixScaleRidgeDetection: ImageFilter 
    {
        #region Fields
        /// <summary>
        /// scale
        /// </summary>
        private double _t;
        private double _epsilon = 10;
        #endregion

        #region Properties
        public double T { get { return _t; } set { _t = value; } }
        public double Epsilon { get { return _epsilon; } set { _epsilon = value; } }
        #endregion

        #region constructor
        /// <summary>
        /// empty constructor
        /// </summary>
        public FixScaleRidgeDetection() : this(1.0) { }
        /// <summary>
        /// imput standard deviation of the distribution
        /// </summary>
        /// <param name="sigma"></param>
        public FixScaleRidgeDetection(double t)
        {
            _flag = "[Image -> Gradient -> Fix Scale Ridge Detection]";
            _t = t;
        }
        #endregion

        #region Apply filter
        //public Bitmap ApplyGraylevelImageFilter(Bitmap srcImage) { return Apply(srcImage, out _log); }
        //public PZMath_matrix ApplyDoubleGraylevelImageFilter(PZMath_matrix srcMatrix) { return ApplyDouble(srcMatrix, out _log); }

        //private Bitmap Apply(Bitmap srcImage, out string logText)
        //{
        //    PZMath_matrix srcMatrix = _converter.ConvertGraylevelBitmapToMatrix(srcImage);
        //    PZMath_matrix dstMatrix = ApplyDouble(srcMatrix, out logText);
        //    Bitmap dstImage = _converter.ConvertMatrixToGraylevelBitmap(dstMatrix, true);
        //    return dstImage;
        //}   

        protected override PZMath_matrix ApplyDoubleGraylevelImageFilter(PZMath_matrix srcMatrix, out string logText)
        {
            // 1. Gaussian Filter
            GaussianSmoothing gaussianSmoothing = new GaussianSmoothing(_t);
            PZMath_matrix L = gaussianSmoothing.ApplyDoubleGraylevelImageFilter(srcMatrix);
            
            // 2. Lx
            GradientXCentralDifference gradientXFilter = new GradientXCentralDifference();
            PZMath_matrix Lx = gradientXFilter.ApplyDoubleGraylevelImageFilter(L);

            // 3. Ly
            GradientYCentralDifference gradientYFilter = new GradientYCentralDifference();
            PZMath_matrix Ly = gradientYFilter.ApplyDoubleGraylevelImageFilter(L);

            // 4. Lxx, Lxy, Lyy
            PZMath_matrix Lxx = gradientXFilter.ApplyDoubleGraylevelImageFilter(Lx);
            PZMath_matrix Lxy = gradientYFilter.ApplyDoubleGraylevelImageFilter(Lx);
            PZMath_matrix Lyy = gradientYFilter.ApplyDoubleGraylevelImageFilter(Ly);

            // 5. Lv, Luu, Luv, Lvv
            int width = L.ColumnCount;
            int height = L.RowCount;
            PZMath_matrix Lv = new PZMath_matrix(height, width);
            PZMath_matrix Luu = new PZMath_matrix(height, width);
            PZMath_matrix Luv = new PZMath_matrix(height, width);
            PZMath_matrix Lvv = new PZMath_matrix(height, width);
            PZMath_matrix dstMatrix = new PZMath_matrix(height, width);

            PZMath_matrix dstMatrixtemp = new PZMath_matrix(height, width);

            dstMatrix.Setall(255);
            //double cosAlpha, sinAlpha;
            double lx, ly, lv, lxx, lxy, lyy, luu, luv, lvv, lv2;
            for (int x = 0; x < width; x++)
            {
                for (int y = 0; y < height; y++)
                {
                    lx = Lx[y, x];
                    ly = Ly[y, x];
                    //cosAlpha = lx / System.Math.Sqrt(lx * lx + ly * ly);
                    //sinAlpha = ly / System.Math.Sqrt(lx * lx + ly * ly);
                    //lv = cosAlpha * lx + sinAlpha * ly;
                    lv = System.Math.Sqrt(lx * lx + ly * ly);
                    Lv[y, x] = lv;
                    lv2 = lv * lv;
                    lxx = Lxx[y, x];
                    lxy = Lxy[y, x];
                    lyy = Lyy[y, x];
                    luu = (lx * lx * lyy - 2 * lx * ly * lxy + ly * ly * lxx) / lv2;
                    luv = (lx * ly * (lxx - lyy) - (lx * lx - ly * ly) * lxy) / lv2;
                    lvv = (lx * lx * lxx + 2 * lx * ly * lxy + ly * ly * lyy) / lv2;
                    Luu[y, x] = luu;
                    Luv[y, x] = luv;
                    Lvv[y, x] = lvv;

                    // detect ridge/valley
                    if (System.Math.Abs(luv) < _epsilon
                        && luu * luu - lvv * lvv >= 0
                        && luu < 0)
                    {
                        dstMatrix[y, x] = 0.0;   // label ridge
                    }

                    dstMatrixtemp[y, x] = System.Math.Abs(luv);
                }
            }
                       
            logText = "succeeds.";
            dstMatrixtemp.Normalize(0, 255);
            //return dstMatrixtemp;
            return dstMatrix;
        } // ApplyDoubleGraylevelImageFilter()

        protected override PZMath_matrix[] ApplyDoubleColourImageFilter(PZMath_matrix[] srcMatrixArray, out string logText)
        {
            // TODO: add this function

            throw new Exception("The method or operation is not implemented.");
        } // ApplyDoubleColourImageFilter()
        #endregion
    }
}
